Nonequilibrium Cooper pairing in the nonadiabatic Regime 
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We obtain a complete solution for the mean-field dynamics of the BCS paired state with a large, 
but finite number of Cooper pairs in the nonadiabatic regime. We show that the problem reduces 
to a classical integrable Hamiltonian system and derive a complete set of its integrals of motion. 
The condensate exhibits irregular multi-frequency oscillations ergodically exploring the part of the 
phase-space allowed by the conservation laws. In the thermodynamic limit however the system can 
asymptotically reach a steady state. 



The study of the dynamics of the BCS superconduc- 
tors has a long history 0. Early attempts to describe 
nonstationary superconductivity were based on the time- 
dependent Ginzburg-Landau (TDGL) equation 0,0,01' 
which reduces the problem to the time evolution of a 
single collective order parameter A(t). The TDGL ap- 
proach is valid only provided the system quickly reaches 
an equilibrium with the instantaneous value of A(i), i.e. 
a local equilibrium is established faster than the time 
scale of the order parameter variation, ta — 1/A. This 
requirement limits the applicability of the TDGL to spe- 
cial situations where pair breaking dominates, e.g. due to 
a large concentration of magnetic impurities. An alter- 
native to TDGL is the Boltzmann kinetic equation^, [|| 
for the quasiparticle distribution function coupled to a 
self-consistent equation for A(f). This approach is justi- 
fied only when external parameters change slowly on the 
ta time scale, so that the system can be characterized 
by a quasiparticle distribution. 

Is it possible to describe theoretically the dynamics of 
a BCS paired state in the nonadiabatic regime when ex- 
ternal parameters change substantially on the ta time 
scale? In particular, an important question is whether, 
following a sudden perturbation, the condensate reaches 
a steady state on a ta time scale or on a much longer 
quasiparticle energy relaxation time scale T e . In the 
nonadiabatic regime both TDGL and the Boltzmann ki- 
netic equations fail and one has to deal with the coupled 
coherent d yna mics of individual Cooper pairs. Recent 
studies[llE3,ElE3 

of this outstanding problem were 
motivated by experiments on fermionic pairing in cold 
atomic alkali gases 0,0|- The strength of pairing inter- 
actions in these systems can be fine tuned rapidly by a 
magnetic field, making it easier than in metals to access 
the nonadiabatic regime experimentally. 

The main result of the present paper is an explicit gen- 
eral solution for the dynamics of the BCS model, which 
describes a spatially homogenous condensate at times 
t <C T e . We employ the usual BCS mean-field approx- 
imation, which is accurate when the number of Cooper 



pairs is large 0|. It turns out that the mean- field 
BCS dynamics can be formulated as a nonlinear classical 
Hamiltonian problem. We obtain the exact solution for 
all initial conditions and a complete set of integrals of 
motion for the mean-field BCS dynamics. 

In this paper we assume that the number of Cooper 
pairs in the system is arbitrary large, but finite. In this 
case the typical evolution at times t <C t £ is guasz-periodic 
with a large number of incommensurate frequencies. The 
condensate exhibits irregular multi-frequency oscillations 
ergodically exploring the part of the phase-space allowed 
by the conservation laws. The system returns arbitrar- 
ily close to its initial state at irregular time intervals. 
However, the return time diverges in the thermodynamic 
limit for most physical initial conditions, while the solu- 
tion asymptotically reaches a steady state on the ta time 
scale. The system thermalizes on a much larger energy 
relaxation time scale T e [t|. 

The dynamics of the BCS condensate following a sud- 
den change of external parameters has been previo usly 
discussed by a number of authors [H E Q E |3 
UtI Hil fl9| . Most notably, a linear analysis around the 
BCS ground state has been performedjllj, [12J and some 
simple particular solutions for the nonlinear mean-field 
dynamics in the context of superconductivity have been 
reportedflj, [l5|. We discuss below how these results fit 
into the general picture. 

We begin our description of the nonequilibrium Cooper 
pairing in the non-dissipative regime, t <C r e , with the 
BCS modelHSEllil. 
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where tj are single-particle energies. The pairing is be- 
tween time reversed states \j j) and \j |)J23]. Our goal 
is to determine the evolution of a state that was driven 
out of equilibrium at, say, t = 0. 

There are several equivalent ways to derive mean-field 
equations of motion. One can start with the BCS product 
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state, Ylj \ Uj(t) + Vj-(i)ct|ct^J |0), and use Bogoliubov- 

de Gennes equations for the time-dependent amplitudes 
Uj(t) and Vj(t). Alternatively, one can study the evo- 
lution of the normal, Gj(t) = —i{[cj^(t),c^(t)]), and 
anomalous, Fj(t) = —i([cj]-(t),Cji(t)]), Green's functions 
at coinciding times [l^. 

The most convenient for us approach to the BCS mean- 
field dynamics is based on the Anderson pseudospin 
representation 1 1 lj . Within this approach the mean-field 
equations are Hamiltonian equations of motion for a 
classical spin chain. Pseudospin-1/2 operators are re- 
lated to fermion creation and annihilation operators via 
K% = {n n + n jl - l)/2 and KJ = c n c 3] 



Pseudospins are defined on empty and doubly occupied 
(unblocked) single-particle orbitals tj. Singly occupied 
orbitals are decoupled from the dynamics. For n un- 
blocked orbitals the Hamiltonian has the form 



n-l 

Hbcs = 2e A Z 
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The mean field approximation is accurate|lCt 111) in 
the thermodynamic limit due to the infinite range of in- 
teractions between spins in the Hamiltonian (J2J. There- 
fore, the effective field seen by each pseudospin in (J2J 
can be replaced with its quantum mechanical average, 
bj(t) = (-2A a .(t),-2A„(t),2c i ), where A(t) ee A x (t) - 
iA y (t) ee gJ2j{Kj '(*)) is the BCS gap function. In this 
approximation, each spin evolves in the self-consistent 

hj x Kj. Taking the quan- 



field: K, 
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turn mechanical average of these equations with respect 
to the time-dependent state of the system, we obtain for 
Bj (t) = (K,(t)) 



n-l 

Sj = bj x Sj bj = (-2gJ x , -2gJ y , 2e 3 ) J = ^ s q 

q=0 

(3) 

The components of the classical spins s| (t) and = s J ± 
is v - are related to Bogoliubov amplitudes and equal times 



Green's functions as 2f 
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Gj(t) — isj[t), Fj(t) — isj(t), respectively. Evolution 
equations yj) conserve the square of the average for each 
spin: dSj/dt = 0. If the spins initially were in a product 
state, s 2 = 1/4. Note also that A(t) = gJ-(t). 

One can check that Eqs. J2J are equations of motion 
for a classical spin Hamiltonian 
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It means that Eqs. J3J are Hamilton equations Sj — 
{Hbcs, s j} derived from Hamiltonian using the usual 
angular momentum Poisson brackets 



(5) 



where a, b, and c stand for spatial indexes x, y, and z. 
The classical model Q can be obtained from its quan- 
tum counterpart (j2J by replacing operators with classi- 
cal dynamical variables and commutators with Poisson 
brackets. 

Both the classical and quantum models and (J2J 
are integrable [23.l25ll2q |. To show this, one can introduce 
a vector function (Lax-vector) of an auxiliary parameter 
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where z is a unit vector along the z axis. Poisson brackets 
between components of L(it) at different values of u can 
be evaluated using Eq. (J5J. 

{L«{v),L\ W )}=e a J C ^- LC{w) (7) 
L ' v — w 

(Relations Q hold for each term in © separately; all 
terms Poisson-commute with each other). It follows from 
Eq. (JZJ) that the lengths of the Lax vector at different 
values of u Poisson-commute: 

{L 2 H,L 2 M}=0 (8) 

The scalar function L 2 (u) can be represented in the form 



(9) 




where 



k=0 
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Since Eq. (JSJ) holds for any v and w, all Hj Poisson- 



commute with each other. Therefore, 



each Hj, as well 



as any algebraic combination of Hj, defines a classi- 
cal model|27j that has n degrees of freedom (n classical 
spins) and n integrals of motion (including itself) and 
thus is integrable in the usual sense [3]. Note that the 
sum of Hj is proportional to the z-component of the to- 
tal spin J, therefore J z is conserved by all Hj and their 
combinations. Moreover, the following identity follows 
from Eqs. (|4TTU|) 



Hbcs = -gJ2 e J H J 



const 



(11) 



This implies that the classical BCS model @ Poisson- 
commutes with all -ffj 's and thus is also integrable. Equa- 
tions l|10(l and (jl 1|) can be straightforwardly quantized by 
replacing Sj — > Kj . The resulting operators Hj all pair- 
wise commute, thus showing the integrability of quantum 
models (QJ and @. 
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To obtain the general solution for the mean-field dy- 
namics of the BCS model Q, we follow the method of 
Ref.Hand introduce n - Reparation variables Uk as ze- 
ros of = L x (u)-iL y (u), i.e. £\ s j/( u k zjy) = 0. 

Equations of motion for the variables Uk are|30| 

ilk = ZiVQfrJuk) Y], ( Uk ~ u^y 1 ( 12 ) 

where Q2n{v) is the spectral polynomial defined as 

Q 2n (u)=g 2 L 2 (u)l[(u-e 1 ) 2 (13) 

j 

By Eq. ©, the coefficients of Q2n(u) depend only on the 
integrals of motion Hj. 

Eqs. (|12|) constitute the well-known Jacobi's inver- 
sion problem solvable in terms of hyperelliptic theta 
functions [3 1|. Here we outline the final answer, the de- 
tails will be reported elsewhere [3f)j. Klenian a- and £- 
functions of genus G (in our case G = n — 1) are defined 
as 

dlno(x) ^ exp[S m (x)/2] 

fe < m ^ ( 14 ) 

<Sm(aO = x • i]uj 1 x + 2z7r(m • rm + to 1 x • m) 

where the sum is over all G-dimensional integer vectors 
m, t — a/a; -1 , and oj, u>' , and 77 are G x G matrices of 
periods (see below). The solution is 

sj(t) = (ZjlWiiit)) = J-(t)r{e 0l t) TJ — ^— (15) 

A(t) = S J_(i) =. ? E S 7W = c - e ^ 4 ^ ± ^ ( 16 ) 

' — ' J (Tlx — a 

i 7 

Here x T = i(c\, . . . , c„_2, 2t + c n _i); d is a vector of 
constants; /3 = gJ z e j! ci, . . . , c n are constants fixed 
by the initial conditions, and 

n-l 

r(u, t) = 1 - 5Z(Cfc(x + d) - C fe (x - d) + a fe ) U fe -" (17) 
fc=i 

Constants a^, the matrices of periods, constant vector d, 
are all uniquely determined[32( by the spectral polyno- 
mial Q2n(u), i.e. by integrals of motion. 

The evolution of Sj(t) = (Kj(t)) described by the gen- 
eral solution is typical of an integrable system [28j|. It 
is characterize d by n frequencies, which in our case can 
be determined|30( in terms of integrals of motion, and 
are typically incommensurate. Note that |A(t)| contains 
only n — 1 frequencies. The typical dynamics is stable 
against perturbations destroying integrabilitvfisj. 

Now let us discuss some particular solutions. There are 
two types of equilibrium states that play an important 



role in the dynamics. In normal states all spins are par- 
allel to the z-axis, 2s| = ±1. Since 2s j = [hj) — 1, these 
states correspond to the ground state and excitations of 
the single-particle part of the Hamiltonian (Fermi 
gas). They are stationary within the mean-field dynam- 
ics |JS). For a finite system, they are non-stationary 
for the quantum Hamiltonian and their short time 
dynamics is entirely driven by quantum corrections (cf. 
Refs.miH- 

The second type of equilibrium states are anomalous 
ones, which correspond to the BCS ground state and 
excitations. These states are obtained by aligning each 
spin in {3J self-consistently along the effective magnetic 
field acting on it. The self-consistency condition is the 
BCS gap equation. As Sj = (Kj), one can obtain the 
BCS wave-function and energy spectrum from anoma- 
lous equilibrium configurations of classical spins s_,-. 

It turns out that equilibrium states are a part of a 
more general scheme when the dynamics of n spins de- 
generates to that of m < 71 collective spins (m-spin solu- 
tions) governed by the same Hamiltonian J3J only with m 
spins and new parameters /ij instead of 6j. Normal and 
anomalous states correspond to 0- and 1-spin solutions, 
respectively. To construct m-spin solutions one has to 
take the Lax-vector (jSJ to be proportional to that of a 
system with m spins tfc, L(u) = + bj/(u — 6j)]Lt(u), 
where bj are time- independent constants, and Lt(ii) 

7, <j + tfc/(u — /ifc). Then, 2(n — m) of 2n, typically 
distinct, roots of the spectral polynomial Q2n(u) become 
doubly degenerate and n — m separation variables Uk are 
frozen in these double roots, which automatically solves 
the equations of motion for these n — m variables. The 
dynamics is obtained by replacing n — > m and e/. — > /j^ 
in Eqs. I|14H17(I and is characterized by m < n typically 
incommensurate frequencies. For m = 2 the solution is in 
terms of hyperelliptic functions of genus G = m — 1 = 1, 
i.e. in terms of ordinary elliptic functions. 

Now let us discuss the connection of our results with 
the previous work. The solutions for the mean field BCS 
dynamics obtained in Ref. are 2-spin solutions in the 
above classification. They were used in Ref. to de- 
scribe the evolution beginning from a state infinitesimally 
close to the normal ground state. In our view, the dy- 
namics in the vicinity of this state can have additional 
features and deserves further analysis. 

The 2-spin solutions resemble the TDGL approach in 
that they describe the dynamics of all pairs in terms of 
only two collective degrees of freedom resulting in large 
amplitude single frequency (periodic) oscillations of the 
order parameter magnitude |A(t)|. Mathematically, they 
lie on a Id curve of points in a multi-dimensional (infinite- 
dimensional in the thermodynamic limit) space of pos- 
sible values of integrals of motion. The situation with 
other few spin solutions is similar pjjj. In contrast, the 
general solution we obtained here typically has a large 
(infinite in the thermodynamic limit) number of incom- 
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mensurate frequencies and a substantially reduced am- 
plitude. The difference between the general and few spin 
solutions is clear in a linear analysis [ill LL=f around the 
BCS ground state that displays normal modes with fre- 
quencies ujk = 2^/e| + Aq, where e& are single-particle 
energies and Ao is the equilibrium order parameter. In 
the linear regime, the general solution becomes an arbi- 
trary superposition of all normal modes, while few spin 
solutions single out all, but few modes. For example, 2- 
spin solutions of Ref. |l5j correspond to a single normal 
mode with a frequency 2Ao. 

In conclusion, we have obtained the explicit general 
solution for the mean-field dynamics of the BCS paired 
state and discussed a number of special cases including 
two types of equilibrium states and few spin solutions. 
A still open problem is to fully analyze the solution in 
the thermodynamic limit. It is also desirable to better 
understand the dynamics in the vicinity of normal states 
where quantum effects become important. Finally, it is 
interesting to identify experimental setups where pecu- 
liar features of the nonequilibrium Cooper pairing in the 
nonadiabatic regime can be observed such as e.g. cold 
Fermi gases. 

We are grateful to I. Aleiner, V. Falko, L. Glazman, L. 
Levitov, A. Millis, A. Polyakov, and O. Tsyplyatyev for 
stimulating discussions. This research was supported by 
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Note Added. Recently, we became aware of a 
publication 36] that is in agreement with some of our con- 
clusions - that the initial dynamics of the normal ground 
state is driven by quantum corrections and that the sys- 
tem can reach a steady state at large times. 
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